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ABSTRACT 



Context. Relativistic jets emerging from active galactic nuclei (AGN) cores transfer energy from the core of the AGN to 
their surrounding interstellar/intergalactic medium through shock-related and hydrodynamic instability mechanisms. 
Because jets are observed to have finite opening angles, one needs to quantify the role of conical versus cylindrical jet 
propagation in this energy transfer. 

Aims. We adopt parameters representative for Faranoff-Riley class II AGN jets with finite opening angles. We study 
how such an opening angle affects the overall dynamics of the jet and its interaction with its surrounding medium and 
therefore how it influences the energy transfer between the AGN and the external medium. We also point out how the 
characteristics of this external medium, such as its density profile, play a role in the dynamics. 

Methods. This study exploits our parallel adaptive mesh refinement code MPI-AMRVAC with its special relativistic 
hydrodynamic model, incorporating an equation of state with varying effective polytropic index. We initially studied 
mildly underdense jets up to opening angles of 10 degrees, at Lorentz factors of about 10, inspired by input parameters 
derived from observations. Instantaneous quantifications of the various interstellar medium (ISM) volumes affected by 
jet injection and their energy content allows one to quantify the role of mixing versus shock-heated cocoon regions over 
the simulated time intervals. 

Results. We show that a wider opening angle jet results in a faster deceleration of the jet and leads to a wider radial 
expansion zone dominated by Kelvin-Helmholtz and Rayleigh- Taylor instabilities. The energy transfer mainly occurs 
in the shocked ISM region by both the frontal bow shock and cocoon-traversing shock waves, in a roughly 3 to 1 ratio 
to the energy transfer of the mixing zone, for a 5 degree opening angle jet. The formation of knots along the jet may 
be related to X-ray emission blobs known from observations. A rarefaction wave induces a dynamically formed layered 
structure of the jet beam. 

Conclusions. Finite opening angle jets can efficiently transfer significant fractions (25 % up to 70 %) of their injected 
energy over a growing region of shocked ISM matter. The role of the ISM stratification is prominent for determining 
the overall volume that is affected by relativistic jet injection. While our current 2D simulations give us clear insights 
into the propagation characteristics of finite opening angle, hydrodynamic relativistic jets, we need to expand this work 
to 3D. 
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1. Introduction Many astrophysical relativistic jets have a high Lorentz 

factor, with 7 yarying from 5 to 30 for AGN jets 

Relativistic plasma jets have been exhaustively studied and (|Kellermann et_aL| |20oJ). Therefore, relativistic dynamics 

. £ observed in detail in radio galaxies such as Abell 2052 models are necessary to study them. 

X ' | Her ^ u les A ' or C yS nus A (|Krausej |2005fc |Boettcher etaj] j e ts associated with AGNs have been separated into 

?H \ 12012J). From multi- wavelength observations in specific different classes. Among the observational classification of 

P_ . sources we have understood important aspects of the ways ra dioloud AG N, we find weaksource F anaroff-Riley class I 

in which jets are emitted by activ e galactic nuclei (AGNs) objects (FR-I. iFanaroff fc Rilevl (jl974D ). with a brighter lu- 



and form pronounced radio lobe s (Blandford fc Reedll974r. mmosity near the ori gin of the jet An example is the jet 
|Scheuer|[l973 |Cohen et al.||2007t |Lister et al.|[2009|). in thc 3C 296 object (TLaing et al.ll200l . Powerful sources 



wu-i u i- , j , are known as Faranoff-Rilcy class II radio galaxies (FR-II) 

While observations give access to dynamic parameters ... . . , , ,. , , * , . 

n ,i t r i c ,, and lobe-dominated radio-loud quasars. An example is the 

such as the Lorentz factor 7 or an estimate of the energy . , . ; , , . . i „ . — , n i h^nh T n ,1 

fl r- ,1 ,,, ■ ., 1 * 1 ■ • • let m the 3C216 object leies et al. 1992). In all cases, the 

flux from the overall luminosity, uncertainties remain m •; , , , , . . . , . . , 

,. f . , ., lf f ■ , ., let has been shown to have a relativistic velocity with a 

the properties of the let matter itself, tor instance, its com- J T r , r _ ,, _,_ T . T , r , 

, , , -4-1 j Lorentz factor of up to 5 for 1R-I and a Lorentz factor 

position, density or pressure and the precise topology and » ^ j, ^ | Cnhe , et a1 J 

strength of the magnetic held. Simulation is a powerful tool 11 r— — ■ — ni r.»r.nK «r ; — 7 , — TT . ' 

to make estimates these parameters, and one can study dif- concentrate here on FR-II jet 

ferent scenarios with them, e.g. the composition of the jet 

consisting of mainly electrons and proton s, or electrons and Previous works have shown how a je t can be launched on 

positrons (jKundt fc Gopal- Krishna 198(j). opening magnetic field configurations ([Blandfor d fc Pavnel 
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1982). For relativistic magnetized winds from a central ob- 
ject, additional collimation by disk winds may be needed 
for th eir ultimate collimation, see iBogovalov fc Tsinganosl 
(2005). Those works illustrate how hard it is to ob- 
tain a fully collimated jet with an overall cylindrical ge- 
ometry. Obs ervations clearly show fi nite opening angle s 
of jets, e.g. iPushkarev et all (|201lh ; IZezas et all (|2005h : 
iFrazer fc Kassiml ^OOoT T Flaring models that fit observa- 
tions well were m ade for FR-I jets, as in the work of 
lLaing et al.l (|2006h . In this model, the original angle ob- 
tained from collimation in the source increases due to in- 
ternal characteristics. For FR-II jets, although the beam 
still shows a finite width and opening angle, this conical 
geometry is often discarded by taking a wide r width of the 
beam , or assuming a precession of the source (|Q'Neill et all 
l2005tl . Few mo dels invoke an actua l conical geometry up to 
a few degrees (|Meliani et al.l [2008). The cylindrical injec- 
tion models dominate because FR-II jets are believed to be 
quickly re-collimated by internal shocks. Because we stud- 
ied jet propagation beyond 0.5 parsec away from the central 
object, we ignored general relativistic as well as magnetic ef- 
fects and performed pure hydro studies for varying opening 
angles. We assumed axisymmetry of the jet and therefore 
used 2D simulations. This allowed us to explore paramet- 
rically hydrodynamic relativistic simulations. Nevertheless, 
a more realistic study would imply 3D simulations, which 
we will perform in forthcoming studies. 

The interest in studing relativistic jets associated with 
AGNs has increased recently, since jets feed back energy 
into the interstellar/intergalactic (ISM/IGM) medium, 
which makes them candidate feedback processes at play 
in the larger scale formation of galaxy clusters. For that 
reason we need to quantify the energy transfer efficiencies 
of AGN jets for cylindrical and conical jets. 

Multiple mechanisms are thought to be responsible for 
the energy transfer between the jet and its surroundings. 
Direct energy conversion takes place at the end of the jet 
beam, immediately beyond the Mach disk surface of the jet 
(in a hydro case), where kinetic energy becomes thermal en- 
ergy through efficient relativistic shocks. Hydrodynamic in- 
stabilities taking place at the interface between the shocked 
jet and shocked ISM matter mix the two fluids, allowing for 
energy transfer from the now hot outer part of the beam 
to the surroundings of the jet. Additionally, a bow shock 
forms around the entire jet cocoon, propagating in all di- 
rections and transferring energy to further regions of space. 
We quantify the relevance of each and how a finite opening 
angle at the jet inlet affects these energy transfers. 

Many other dynamical features are known from observa- 
tions, for instance the existence of X-ray blobs along the jet 
path. Different scenarios are studied to explain this emis- 
sion, such as particle acceleration at bow shocks when stars 
cross the beam of the jet. We propose that our jet dynam- 
ics model to explain those blobs without invoking external 
agents such as star-related bow shocks. This explanation 
revisits the role of reco nfinement shocks to explai n these 
knots, as in the work by iKomissarov fc Falld (|1997f ). 

In the second part of this article we describe the rela- 
tivistic hydrodynamic model and the way we use it with 
our code, MPI-AMRVAC. We discuss our settings and ini- 
tial conditions and their motivation. Part 3 will compare 
our different simulations by giving a parametric study of 
the effect of an opening angle on the dynamics of our jets. 
This comparison studies the geometry and energy content 



of the jets. The fourth part discusses two aspects found in 
our simulations that relate to the internal structure of the 
jet: a scenario for X-ray blobs and the dynamically formed, 
layered structure of the jet. We finish with general conclu- 
sions obtained from this work. 

2. Simulation 

2.1. Governing equations 

Like in the previous work of Meliani e t alJ 



we assumed axisymmetry. We use a relativistic hy- 
drodynamic model one with the relativistic variant of 
the Euler equations with th e Mathews approximation 
(jBlumenthal fc Mathewill976l ) to the Synge gas equation 
as a closure. They express the conservation of the different 
relevant quantities. The first equation is the conservation 
of particle number 



dpi 
dt 



^ + V-(p 7 ^) = 0, 



(1) 



where 7 is the Lorentz factor, p is the proper density and v 
is the three- velocity of the fluid in the chosen lab frame. The 
second equation is the conservation of momentum, written 
as 



dS_ 

dt 



V -(Sv) + Vp = 0, 



(2) 



where p is the pressure in the fluid frame and S is the 
momentum density defined as 



S = hj 2 pv . 

The relativistic specific enthalpy h is defined as 



(3) 



(4) 



where T is a parameter, taken as the non-relativistic value of 
the polytropic index of the gas | , e = m p +eth is the specific 
internal energy including the rest (proton) mass m p and 
the specific thermal energy e t h- This is equivalent to using 

a spatially varying, effective r c fr = T — ^p-(l — ^r) that 
varies between its classical value of 5/3 and its relativistic 
value of 4/3. The third equation expresses the conservation 
of energy density 



dr 



V ■ (5 - (rfv) = . 



(5) 



where the energy density r in the lab frame with the rest 
mass contribution already subtracted, is written as 



r = ph-f 2 —p — (r{. (6) 

The Mathews approximation to the Synge gas equation is 
our closure relation, and it writes as 



V 



P\- f" 



(7) 



The speed of light is taken equal to unity in our model and 
therefore in all equations above. 
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2.2. Code and implementation 

The code we used is MPI-AMRVACQ (jKeppens et al.ll2012ft . 
designed to solve a set of equations in conservative form, 
as is the case for our relativistic HD model. The solver is a 
hybrid version of the HLLC scheme, for which a validation 
can be found in the work of lMeliani et al.l (|2008[ ) . The nov- 
elty of this scheme is that it switches to a TVDLF solver 
where the HLLC solver could induce spurious oscillations. 
More validation tests for relativistic HD an d MHD can be 
found on the code website and in the work of lKeppens et all 
(|20T1 . 

The code incorporates block-tree adaptive mesh refine- 
ment (AMR), and we used a weighted combination of the 
normalized second derivative of the density variable p7, the 
pressure p and the Lorentz factor to follow the dynamics of 
the jet and automatically refine the grid where dynamics 
takes place. Additionally, a forced refinement of the grid 
was activated around the inlet region of the jet, up to the 
maximum level specified, to accurately handle the bottom 
boundary. 



Beam luminosity 


Domain size (pc) 


Domain 
resolution 


10 4b erg.s~ L 


[0.5 20]x[0 5] 


2400x1520 


Beam radius (pc) 


Beam 7 




0.05 


10 





Table 1. General parameters for all simulations. 



Case 


Density ratio 


Opening 
angle 


Density 
profile 


A 


0.01 


5° 


uniform 


B 


0.0001 


5° 


uniform 


C 


0.018 


0° 


King 


D 


0.018 


5° 


King 


E 


0.018 


10° 


King 


F 


0.1 


0° 


King II 


G 


0.1 


5° 


King II 


H 


0.1 


10° 


King II 



Table 2. Parameters distinguish between the simulations. 



2.3. Initial conditions 

The rectangular domain was defined to be 20 by 5 pc jet 
beam radius in length and width, respectively (z, r), with a 
length normalization of 1 parsec which corresponds to 400 
by 100 jet beam radius as we show below. The domain did 
not include the central engine that launches the jet, but 
started at z = 0.5 distance from the source, therefore we 
can avoid treating the close surroundings of an AGN core, 
which requires general relativistic theory. From a numerical 
point of view, all AMR sets start with a base resolution of 
300 by 190 grid points, in (z, r) respectively, allowing four 
levels of refinement (hence an effective resolution of 2400 
by 1520). We made sure that the resolution allowed us to 
resolve the jet properly: the jet width on inlet, fixed to be 
r& = 0.05 parsec, is 15 grid points wide at the maximal level 
of refinement taht is enforced in the jet region. We set the 
density on the inlet of the jet pj e t such that the jet total 
kinetic luminosity is Lj Ct — 10 46 erg.s -1 with a Lorentz 
factor 7 equal to 10, in a jet with radius rj, given above, all 
values representative for an FR-II type AGN jet. 

For the ISM we chose a model sim ilar to the beta model 
of the King atmosphere (|Kind[l962T ) profile. This assumes 
that the ISM has been structured into a den sity variation 
by the winds comin g from the central object (jO'Neill et al.l 
[20051 1k raus i2005f) . resulting in a slowly decreasing density 
inside a core radius. Outside of this radius the density drops 
rapidly up to an asymptotic value. This ISM density then 
follows 

( z 2 Y' 3 
Pism = Po ( 1 + , (8) 

where z c is setting the size of the core of the galactic source, 
po quantifies the density at the source of the jet (p = 5000 
g.cm -3 for cases A to E and po = 500 g.cm -3 for cases F 
to H), and f3 is a parameter set here to 0.75. The distance 
z is here the axial direction instead of the actual radial 
distance from the central engine. This was done for numer- 
ical convenience, because it more easily handles the bot- 
tom boundary and avoids an artificial 'pixelization' of the 
density background at low refinement levels. We moreover 

1 Code homepage at 'homes.esat.kuleuven.be/~keppens' 



assumed that the jet will not propagate over a large radial 
distance, resulting in a distance to the source d source ~ z 
for r small. Normally, z c is set to a few kpc, corresponding 
to the average size of a galaxy. This would result in a flat 
density profile in our parsec scale simulation. This is the 
case for simulations A and B, which were injected in uni- 
form media, denoted 'uniform' in Table [U For simulations 
C, D and E we set z c to be 0.5 to investigate a possible 
change in the propagation of the jet as it becomes more 
dense than the ISM. With this setting in the domain the 
density profile behaves as a simple decreasing power law. 
This is denoted as 'King' in the summary Tabled Finally, 
simulations F, G and H are in between with a z c at 5 pc, 
called 'King II' in the summary table. 

The non-dimensional numeric value for the pressure was 
set on the whole domain to be uniformly small (po — 10~ 3 po 
for A to E and po = 10 _2 po for cases F to H) so that it 
does not dominate the early stage of the dynamics. Having 
a uniform pressure also removes the need to have a gravity 
potential to keep the ISM in force balance. We call 77 the 
density ratio between the jet and the ISM, i.e. Pj e t/ Pism = V 
at the position z = zq. 

We chose for cases A and B a ratio 7/ equal to 0.01 
and 0.0001, respectively, investigating the dynamics of an 
underdense jet. For cases C, D and E, this ratio was taken 
to be T] = 0.018. For those three cases the density of the 
ISM decreases rapidly in the domain and we expect the 
density of the jet and the ISM to match around the end of 
the domain. Cases F, G and H have a ratio 77 = 0.1 and 
a 'core' inside the domain, resulting in our simulation box 
in a less rapidly decreasing density profile of the ISM, the 
density of both jet and ISM match around the middle of 
the domain. 

We simulated on domains z <E [zq, z out ] and r G [0, r ou t]- 
The boundaries at z out and r out , as seen in figure [TJ are 
open boundaries. The bottom boundary at zq = 0.5 pc is 
an open boundary for r > r\> and was set to the jet injection 
for r < r&, where rt, = 0.05 pc. The r — boundary is the 
symmetry axis. We typically stopped the simulation before 
the jet reached the boundaries z out and r out to prevent 



3 



Remi Monceau-Baroux et al.: Effect of angular 



To Tout 




Central engine 



Fig. 1. Geometry of the jet as used in the simulation: the 
jet enters the domain at zq = 0.5 pc from the source with a 
width setting r\, as calculated assuming a first angle of 5°. 
A flaring-like case assumes a fixed angle 6, which is added 
to a reference cylindrical jet, thereby introducing a 'virtual 
source location'. The actual domain of the simulation is 
comprised between zq and z out for the axial direction and 
r = and r out for the radial direction. The velocity field 
is initiated in a finite region up to z%, according to the 
formulae given in the text. 



energy from leaving the domain because we need to quantify 
energy deposition. 

Fig. [1] gives the geometry of the jet generated with a 
flaring-like scenario. In the initial conditions, the jet is con- 
sidered to already have propagated for sometime and to 
have reached z\ = 1 pc. The initial radius of the jet at inlet 
is fixed at r\> = 0.05 pc and considers the jet to be coming 
out of a central source located at 0.5 parsec below the grid 
starting at zq. Then according to the flaring- like scenario, 
introducing a fixed angle 9, we can obtain the position of 
a virtual source at z v , which for the angles of 5 and 10 de- 
grees is somewhere in between z — and z — Zq, as seen 
in Fig. [TJ This was used to initiate the entire jet region 
within the computational domain, by considering the flux 
coming from this source and adopting a constant mass flux 
along the section of the cone. This means that we initial- 
ized density and velocity in this part of the domain with 
the following equations: 



p(r,Z,t = 0) =Pjet7 —, ^~TT 7^TW> (9) 

(r b + {z- z )tan(6»))^ 



v z (r, z,t = 0) = v cos(arctan( )), (10) 



v T (r, z, t = 0) = v sin(arctan( )). (11) 

z — z v 



opening on the dynamics of relativistic hydro jets 




Fig. 2. Example of our method to instantaneously identify 
jet regions for a jet opening angle of 10° at time t=20 
(case E). We show the masks locating the jet (left) mixing 
region (center) and the shocked ISM (right). We identify 
the jet beam by 7 > 2, finding the mixing region where 
|(V x t))^| > 1 without the jet beam, and the cocoon where 
t =^= Tt—o with beam and mixing region excluded. 



3. Results -parametric study 

3.1. Dynamics - general description 

As previously pointed out, we aim to study the transfer of 
energy between the jet and the ISM. To do this, we need 
to trace and quantify the significance and role of the differ- 
ent regions of the domain that characterize the jet beam, 
the instability-dominated region, the shocked ISM (SISM), 
and the unperturbed ISM. The main feature of the region 
where the beam of the jet is located is its high value for the 
Lorentz factor. We empirically defined the threshold value 
to be 7 > 2. With this criterion, the beam is well-captured 
without pollution from the cocoon region, which still has 
relativistic flows where 2 > 7 > 1. The instability domi- 
nated region was located with the curl of the velocity field: 
because we did not include jet rotation and adopted jet ax- 
isymmetry, only the toroidal (cj>) component of the curl is 
non vanishing. We can instantaneously capture this mix- 
ing region where instabilities develop by locating where the 
magnitude of the curl component is higher than a threshold 
value of 1. The shocked region is identified by quantifying 
where the instantaneous energy density r is not equal to 
its initial value. An example of these criteria applied to the 
simulation of case E at time 20 of the simulation can be seen 
in Fig. [5J Because the regions are not suppose to overlap, 
we substracted the jet from the mixing region, and both jet 
and mixing region from the shocked ISM. 

Figs. [31 and U show instantaneous comparisons between 
cases C, D and E, with a 0°, 5°, and 10°opening angle with 
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rho 




Fig. 3. Comparison of the density distribution at time t=20 
for an opening angle of 0° (left - case C) and 10° (right 
- Case E). A wider opening angle increases the width of 
the jet beam that stays collimated. The maximal reach of 
the jet beam is inversely proportional to the opening angle. 
The radial reach of the shock surrounding the jet remains 
of the same order. 



a King-type density profile, respectively. Common features 
can be described: the beam of the jet first shows an un- 
perturbed region in which the matter propagated without 
interaction. It is then re-collimated and a shock forms that 
finally reaches the axis. Figure O shows the density varia- 
tion along the symmetry axis at the same time t — 20 for 
the three cases C, D and E. For the finite opening angle 
cases D and E, this plot shows that the position of this 
re-collimation shock (at about z w 3) slightly changes for 
varying opening angle. We discuss in section [4. II how these 
shocks, which are repeated along the symmetry direction, 
can provide a scenario for X-ray blobs. 

The head of the jet beam sweeps up matter and the 
wider the opening angle of the jet, the wider this head is and 
the more matter is pushed away. This qualitative statement 
can be quantified using the masks introduced earlier: we can 
estimate the mass swept up by quantifying 

M swept = // pjsu(r, z,t = 0)2nrdrdz 7 (12) 

J J bcam(i) 

which is shown in Fig.[S] A bow shock surrounds the whole 
structure and propagates into the ISM, transferring energy 
to an expanding region. By using the filters previously dis- 
cussed, we performed a series of diagnostics: by integrating 
over the different regions, we can quantify the effect of the 
variation of the opening angle on, e.g. the fraction of the do- 
main volume occupied by the mixing region, which we show 
in Fig. [71 We can also instantaneously locate the overall po- 
sition of the frontal bow shock as seen in the right panel of 
Fig. [5J Therefore we can quantify the maximal radial and 




Fig. 4. Comparison of the proper density distribution at 
time t=20 for an opening angle of 5° (left - case D) and 
10° (right - Case E). Conclusions from visual investigation 
agree qualitatively with fig. [31 



Density cut along the Z axis 




z 



Fig. 5. Density profile along the symmetry axis of the jet 
for cases C, D and E and for the unperturbed atmosphere. 
All cases show an increased density at the head of the jet 
beam. The maximal z reach is inversely proportional to 
the opening angle. Case C (cylindrical case) keeps a near 
constant density. Cases D and E show a series of nodes with 
increased density. The position of the nodes depends on 
the opening angle: the wider the opening angle, the longer 
the internodal distance. The jet density in case C becomes 
higher than the ISM density from about z ks 12, whereas 
for the other two cases, it remains below the ISM density 
profile throughout the domain. 



axial expansion of this region over time as seen in Figs. [S] 
and [5J At the dynamically evolving contact interface be- 
tween matter from the jet and the swept-up ISM, instabili- 
ties develop. These are of two kinds: Rayleigh- Taylor insta- 
bilities forming at the accelerated and decelerated contact 



5 




Fig. 6. Matter swept up by the jet beam for different cases. 
For the same density profile (cases C, D and E, or cases, 
F and H) wider jets sweep away more matter. We can see 
that the density profile clearly matters. 



interface of the jet as it propagates into the ISM, which then 
are advected along the side of the jet. Additionally, Kelvin- 
Hclmholtz instabilities form along the jet beam flanks due 
to the velocity gradient between the moving matter of the 
jet that propagates axially and the slower matter from the 
interaction region with a more complex velocity field. The 
formation of both those instabilities can be followed with 
the filters discussed previously. We quantify how much of 
the energy transfer takes place when these instabilities mix 
matter from the jet and from the ISM. By looking at a plot 
of the curl of the velocity field, we can highlight those insta- 
bilities. As expected, at the interface between the matter 
of the jet with a high Lorentz factor and the more complex 
and slower matter from the interaction region, we can see a 
strong shearing of the fluid with high negative value of the 
curl of the velocity field. We can differentiate likewise the 
unperturbed part of the beam, and the part of the beam 
where recollimating shock waves are found. The front of 
these shock waves, as clearly identified from a pressure plot 
(not shown), coincide with a locally strong shearing. The 
vorticity mapping makes us able to track the instabilities 
and their evolution in time. By integrating the instanta- 
neous energy density over the different regions as defined 
by the filters, we can obtain an estimation of energy transfer 
to the different regions. 

3.2. Dynamics - overall propagation distances 

As mentioned before, Fig. [5] shows the axial offset position 
of the bow shock front over time for an opening angle of 
the jet of 0° (case C), 5° (case D) and 10° (case E) with a 
King- type density profile as defined in section 12751 It shows 
the same for the pair of 0° (case F) and 10° (case H) with a 
King-II-type density profile. Finally, it quantifies the same 
position for an opening angle of 5° with a flat density pro- 
file (case A). All data sets start similarly as the head of 
the jet propagates at early time with the velocity obtained 
fro m the expression for pressure-matched jet propagation in 
ID ([Marti et alj|1997l ). Case A highlights how the density 
structure of the ISM affects the propagation. This uniform 
density case decelerates the head of the jet more as the 
influence of the ISM augments, in contrast to cases C, D 
and E, where the ISM density decreases. The inferred ve- 



Fig. 7. Volume fraction of the mixing region up to a time 
t=20 for different cases. The volume depends both on the 
opening angle and the properties of the ISM density profile. 
The volume increases for a smaller opening angle although 
the radial reach of this region diminishes. The maximal z 
reach of the jet is therefore dominant over the radial ex- 
pansion of the mixing region. 
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5 degree - Flat - Case A 






degree - King - Case C 






5 degree - King - Case D 






10 degree - King - Case E 






degree - Kingll - Case F 






10 degree - Kingll - Case H 
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10 20 30 40 50 60 

time (year) 

Fig. 8. Radial expansion of the shocked region over time for 
different cases. Despite initial differences, the radial width 
of the jet-influenced region does not vary much for a varying 
opening angle. 



locity of the head of the jet, quantified as the slope of the 
curves shown in Fig. [5J has a low Lorentz factor around 
1, and the velocity within the beam of the jet itself has a 
Lorentz factor of 10 or more (due to a compression of the 
flow, which can increase the Lorentz factor above its initial 
value). Different effects influence the acceleration or decel- 
eration of the head: the principal one is the continuous ram 
pressure force exerted by the matter from the jet. As the jet 
beam propagates, the ram pressure force at the Mach disk 
is proportional to its surface area and depends on veloc- 
ity. A significant factor is the accumulated mass of matter 
swept up by the jet beam. This was quantified earlier in 
Fig. [HI We discussed before that one of the main mecha- 
nisms of transfer of energy between the jet and the mixing 
region is dependent on the conversion from kinetic to ther- 
mal energy at the Mach disk, which is then mixed with the 
surroundings of the jet. Therefore, monitoring the velocity 
of the Mach disk, and the closely related propagation speed 
of the frontal bow shock as inferrable from Fig. \Q gives us 
an indication on the efficiency of this energy transfer: one 
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50 60 



Fig. 9. Position of the bow shock marking the head of the 
jet over time for different cases. The density ratio and the 
opening angle are both relevant for the propagation of this 
bow shock. The wider the opening angle, the more it decel- 
erates. The higher the density ratio rj, the slower the shock 
front deceleration. 



Fraction of the domain occupied by the Shocked ISM region 



5 degree - Flat - Case A 
degree - King - Case C 
5 degree - King - Case D 
' " 1 igree - King - Case E 
degree - Kingll - Case F 
1 degree - Kingll - Case H 




Fig. 10. Fraction of the domain volume occupied by the 
shocked ISM region up to a time t=20 for different cases. 
This volume fraction of the shocked region both depends 
on the opening angle and on the density profile of the ISM. 
We know from Fig. [5] that all simulations have comparable 
radial reach, therefore we expect the volume of the cases to 
be ordered as in Fig. [5J This is the case within the same 
series of simulations (cases C, D and E, and cases F and 
H). 



would expect that the kinetic energy from the jet would be 
transferred into thermal energy for the cocoon. 

The upper limit of propagation is a ballistic behavior 
of no deceleration of the jet. The closest to this is the fully 
cylindrical jet of case C (0°, King), consistent with the least 
swept-up matter shown in Fig. [SJ The ram pressure felt 
by the head of the jet as it propagates decreases as the 
density of the ISM decreases in the axial direction. For cases 
D and E (5°and 10°King), we see a clear deceleration of 
the velocity of the head of the jet, due to a wider opening 
angle of the jet. Ultimately, we expect the mass of material 
piled up at the head of the jet to be high enough to stop 
its propagation, and this could happen sooner for wider 
opening angles and lead to shorter lengths of propagation 
(but this is on time-length scales beyond our simulations 
here). The test case scenario A with a flat density profile 



shows a strong deceleration as the jet encounters a constant 
high ISM density. As for the previous cases, the total mass 
of matter pushed aside by the jet increases, decelerating 
the head of the jet even more. This can be seen in Fig. [5] 
where case A (5°, uniform density profile) shows the lowest 
curve. By comparison with case D (still 5°but King density 
profile) it becomes evident that the density profile of the 
ISM is a governing effect. All these results are consistent 
with Fig. [5] a wider opening angle increased the mass of 
matter pushed away and therefore implies a higher value of 
the ram pressure. For flat density profiles, the density does 
not decrease axially and we reach even higher values for 
the swept-up matter and therefore a much more significant 
deceleration of the head of the jet. From Fig. [5] we can 
expect the same relation between the opening angle for the 
jet and the increase of energy transfer between the jet and 
the ISM. We also infer that a wider opening angle for the 
jet translates into a reduced reach of the jet into the ISM. 
This calls for longer simulations to study the longer term, 
far distance behavior. 

Cases F and H (0°and 10°, King II density profile) are 
given in all figures as the limit cases of the second series 
of simulations entering a varying density atmosphere. The 
intermediate opening angle case G (5°, King II) is not quan- 
tified in the figures, but we verified that it is always inter- 
mediate between cases F and H. The conclusions obtained 
from the study of the second serie and their comparison 
between one another are similar to the ones obtained from 
the related series of cases C, D and E (0°, 5°, 10°, King). 
Case B (5°, uniform density profile, r\ = ICP 4 ) is found to 
follow all trends that were identified with case A (5°, uni- 
form density profile, r\ = 10 -2 ), only shifted in intensity 
due to the difference of the density ratio r\. 

In summary: A wider opening angle results in shorter 
axial propagation, while the density of the medium in which 
the jet propagated plays the expected role: a higher density 
of the external medium results in a higher deceleration. 

3.3. Dynamics study - volumes and radial expansion 

As the jet propagates into the ISM, a bow shock front forms 
around it and propagates. This front stays always fairly 
close to the end of the jet beam, but radially it expands 
much faster than the jet which is re-collimated. As this 
bow shock propagates into the ISM, it transmits energy to 
the local matter. This is the main way to transmit energy 
to far radial distances. Fig. [S] shows time evolution of the 
maximal radial distance reached by the SISM bounded by 
the bow shock, for several representative cases. 

For propagation into a King atmosphere, we find that 
this maximal radial expansion of the SISM stays on the 
same order irrespective of the initial jet opening angle. 
While the jet beam on input is only 0.05pc for all cases, the 
bow shock's sideway influence region increases to 5 parsec 
within the simulation time. For an opening angle of 5°, but 
with a flat density profile instead of a King-type atmosphere 
(case A), the radial expansion is slower and stays within 2 
pc in the same timespan. We can easily explain this behav- 
ior by considering that this shock is propagating sideways 
with a speed influenced by the local sound speed. As said 
in section |2~51 the initial pressure in the ISM is constant, 
whereas for case A, density is also constant, but density 
decreases with axial distance for other cases. Therefore the 
local sound speed (which is a factor of ^) in the ISM in- 
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creases in the axial direction. This in combination makes 
the sideways expansion more dependent on the ISM cir- 
cumstances than on internal beam properties. The latter, 
as explained previously, does affect the beam axial propa- 
gation behavior, as shown in Fig. [51 where smaller opening 
angles meant more deceleration. Associated with this dif- 
ference in reach in the axial direction for the jet, Fig. [TU] 
shows a higher volume fraction for the shocked region for 
smaller opening angles for the same density profile (com- 
pare in particular cases C, D and E, i.e. 0°, 5°and 10°, and 
King- type). This figure quantifies the ratio 



Vsis 



M 



fit 



SISM(t) 



2tt r dr dz 



//do 



2n r dr dz 



(13) 



In reality, the pressure in the ISM is not constant as we 
move away from the source. Therefore to link those state- 
ments to observations, we would need to study the evolution 
of the ratio jj of the surroundings of AGN jets, which is still 
a problem for future study. 

From Figs. [3] and UJ one would expect the volume of the 
mixing region for cases C, D and E to be more significant 
for wider opening angle. But from Fig. [JJ which quantifies 
the instability dominated volume fraction, we can see that 
the radial expansion of the instability-influenced region is 
compensated for by the lesser axial propagation of the jet. 
By comparison, case A, which encounters a much higher ex- 
ternal density, and therefore a much higher ram pressure, is 
confined to a smaller volume. For the jet beam region itself 
(as identified at each time with the mask as in the left panel 
of Fig. [5]), the main factor for the variation of its volume 
(figure not shown) is the radial expansion set by the initial 
opening angle: for an opening angle of 0° (case C), the jet 
does not expand radially and therefore its volume increases 
linearly with its axial expansion. For a wider opening angle 
of 5° and even more for 10°, the jet beam has a varying ra- 
dial expansion due to rccollimation by internal shocks and 
the varying ram pressure at its head. 

In conclusion: The radial propagation of the frontal bow 
shock is dominated by the external medium parameters. 
This shock is the main means of energy transfer to large 
radial distances. 

3.4. Temporal evolution of the energy content 

We now discuss the evolution of the energy content by us- 
ing various masks as described earlier in the text, to dif- 
ferentiate the contributions of different spatial regions (see 
section 13.11) . We will do so by quantifying the fraction of 
excess energy present in each part of the domain as 



E, 



mask 



//; 



mask(i) 



3 (r, z, t) 2irr dr dz 



E, 



domain 



//do 



3 (r, z, t) 2irr dr dz ' 



(14) 



where as excess energy we define the energy in each point 
of the domain, minus the contribution from the ISM at the 
same point for t = 0: 



3 (r, z, t) = r(r, z, t) - n SM (r, z,t = 0). 



(15) 



For various choices of masks, we give the temporal evolution 
of the energy fraction present in the different parts of the 
domain, namely the mixing region (Fig. Ill)) , the shocked 
ISM (Fig. [T2D and their combination (Fig. USD. Note that 



the denominator is equivalent to the time integration of 
the jet influx after subtracting a small loss term -Ei oss that 
leaves through the open part of the bottom boundary: 



-^domain — / -^jct dt E\ oas 



3 (r, z, t) 2irr dr dz . 



(16) 



domain 



It is in particular non-zero at t = due to the finite region 
of the domain where we initialize the jet properties (z £ 
[zq,zi] and r within the specific conical segment). For the 
masks mentioned (mixing region and SISM) , the numerator 
starts at a zero value. 

For Figs. QT]- HH the initial behavior up to time 1 in code 
units is mostly due to the way we initiate the simulation 
with a sharp-edged conical jet segment within the domain. 
More relevant for our discussion is the trend of these curves 
after this time. 

In section 13.21 we postulated that the energy transfer 
in the mixing region is dominated by the energy that is 
transferred through the end of the jet beam (Mach disk 
location). Therefore, a relation should exist between the 
total mass of swept up matter by the beam (as quantified 
earlier in Fig. [5]) and the energy transferred to the mixing 
region. Indeed, as we can confirm from Fig. 111! within the 
same series of simulations (cases C, D and E: 0°, 5°and 10°, 
King or F and H: 0°, 10°, King II), by keeping the density 
profile fixed, the wider the opening angle, the higher the 
mass of swept-up matter and the excess energy transfer in 
the mixing region. In all cases, we observe a decrease of the 
energy transfer rate to the mixing zone up to an asymptotic 
regime where the transfer becomes time-independent, e.g. 
cases E and H. Note that there is a clear trend for higher 
energy transfer for larger opening angles, minimizing the 
role of the ISM stratification. 

Fig. [T2] shows the same quantification as before, with 
the mask now identifying the SISM regions. Similarly as 
before, wider opening angles result in higher energy trans- 
fer. In contrast to the mixing region case, the ISM stratifi- 
cation plays an important role in the energy transferred to 
the SIMS region. Indeed, Fig. [T2] shows a clear difference 
between cases C, D and E, and cases F and H. We point 
out that the energy transfer in the shocked ISM region is 
not dominated by the volumetric evolution of this region 
as seen in Fig. [TDJ We can conclude that the dominant pa- 
rameter for energy transfer is here the density profile of the 
ISM. 

Combining both regions, Fig. [T3] then quantifies how 
much energy is transferred to the jet beam surroundings 
as a whole. Contrary to our initial expectation that most 
of the energy transfer should occur at the Mach disk and 
mixing region alone, Fig. (T3] shows that the total energy 
transfer is dominated by the transfer in the SISM region. 
We observe, on average, a ratio of 3 to 1 between the energy 
fraction present in the SIMS region, and its mixing region 
counterpart. Moreover, we identify the obtained percentage 
of the transferred energy as an estimate of the time-average 
energy transfer rate for the duration of our simulations, 



transfer 



(t)dt E, 



transfer 



(*) 



Edo 



At) 



x L 



jet • 



(17) 



In case D (5 degree opening angle with a King density 
profile) we observe a transfer rate of 2.5 x 10 45 erg s _1 (i.e. 
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1 degree - King - Case E 
degree - Kinglh-Case-F — 
10 degree -Kingl! - Case H 




Fig. 11. Fraction of the total domain energy excess present 
in the mixing region, which is dominated by instabilities, 
for the different cases. 



Energy fraction in shocked ISM region 



5 degree - Flat - Case A 
degree - King - Case C 
5 degree - King - Case D 
1 degree - King - Case E 
degree - Kingll - Case F 
1 degree - Kingll - Case H 




Fig. 12. Excess energy ratio in the shocked fSM region for 
different cases. 



Fig. 13. Total excess energy ratio in the mixing region 
and shocked ISM for different cases. By comparison with 
Figs. [TTI and [T2l we conclude that the transfer is dominated 
by heating from the shocks in the shocked ISM region. 



Thermal energy fraction 



5 degree - Flat - Case A 
degree - King - Case C 
5 degree - King - Case D 
1 degree - King - Case E 
degree - Kingll - Case F 
10 degree - Kingll - Case H 



Fig. 14. Thermal term of the total energy transfer. 



at t = 20, 25 % of the total kinetic energy flux -©domain 
is injected into the surrounding medium). We can quantify 
the ratio over thermal and kinetic energy to the total energy 
even better by noting that the local energy density minus 
rest mass, r, can be approximately expressed by its thermal 
and kinetic contribution as 



is dominated by the heating by shocks which happens in the 
SISM region. We expect the transfer by mixing of the jet 
material and the surroundings to be dominant only at later 
times. It is also possible to quantify an average transfer rate 
of energy from the jet to its surroundings. 



^thermal = p(-f 2 „ ° : ~ 1), (18) 
1 off - 1 

Tkinotic = 7(7 - !)PC 2 7 (19) 

where 

r - 1 ml 

r cff = r-— - (i--£) (20) 

is the local varying effective polytropic index of the gas. In 
Fig. 1 141 we use these individual contributions in the integral 
expressions, subtracting their initial ISM values. We then 
find that about 80 percent of the energy transferred to the 
ISM is composed of the thermal term, which again con- 
firms the dominant contribution of the bow-shock-related 
heating. 

In conclusion: it is possible with our code to track where 
the energy is deposited and in which form, kinetic or ther- 
mal. For early stages of the propagation, the energy transfer 



4. Dynamical details in conical relativistic HD jets 

4.1, Internal structure of the jet 

In Fig.0 we compared the density profile along the symme- 
try axis for several cases. By overplotting the ISM density 
profile, it is possible to observe the evolution of the density 
ratio T] = yOjet/ptSM- We recall that this ratio is set to 0.018 
at the bottom boundary for several cases shown. For an 
opening angle of zero, the jet is cylindrical from the start 
of the domain and its density stays constant on average, 
whereas the density of the ISM decreases. Small fluctua- 
tions due to different internal shocks can be observed. For 
this case C (0°, King), after t=10 the density of the jet be- 
comes higher than that from the ISM at a distance of 11 
parsec from the central engine. 

When finite opening angles are present, a clear recolli- 
mation shock forms. Up to this recollimation shock, for an 
opening angle of 5° and 10° the jet density decreases as a 
factor of a , 1 ,„ yi with z the distance to the source and 9 
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the fixed opening angle of the jet. As z increases, the rec- 
ollimation shock wave propagates into the beam and com- 
presses its density. The position at which this shock first 
reaches the symmetry axis evolves in time, until it finds an 
equilibrium recollimation position at 3 parsec for an open- 
ing angle of 5°, and 4 parsec for an opening angle of 10°. 
From time t=10 this quasi-equilibrium position is reached 
for both cases D and E, as in Fig. (similar trends emerge 
for the series F, G and H). A two-dimensional view of the 
near-stationary nodal structures is shown in Fig.[T5]for case 
E at time t = 20. This panel shows the density and effec- 
tive polytropic index near the injection zone. Even if this 
does not dominate the dynamics, magnetic fields should be 
present in reality, enabling processes such as Fermi-11 accel- 
eration of particles in the nodes along the jet. This would 
result in an increase of the radiation processes at those po- 
sitions, which is a scenario to explain the X-ray blobs found 
in observations of AGN jets. 

We have a Synge-like closure relation, defined such that 
the effective polytropic index can vary between 5/3, the 
classical polytropic index for a proton gas, to 4/3, its rela- 
tivistic value. We can recover the local value of this effective 
index and plot it following the equation [20] This allows us 
insight into the thermodynamic condition of the gas. Figs. 
[TBI and [T71 show this effective polytropic index for cases C, 
D and E. We can see for both cases with a finite opening 
angle that the gas expanded and becomes cooler. This is 
due to adiabatic cooling in a freely expanding gas. The re- 
collimation shock reheats the gas and moves its polytropic 
index to the relativistic value of 4/3. A zoom-in of this re- 
gion was shown in Fig. [T5] for case E. The main jet beam 
remains at a relativistic temperature. For an opening an- 
gle of 0° , the jet matter does not experience the adiabatic 
cooling phase and stays uniformly thermodynamically rel- 
ativistic. We note that the core of the instability-related 
vortices also coincide with a mildly relativistic value of the 
polytropic index. In the variation of the effective polytropic 
index, we can also see how several secondary shocks traverse 
the SISM cocoon region, adding to the overall reheating. 

In conclusion: using of an opening angle for the injection 
of the jet modifies its internal structure and gives rise to 
recollimation nodes along its path. 

4.2. Dynamic formation of a layered structure 

We also point out another detail of the finite opening angle 
cases by showing aanother zoomed view on the near inlet 
expansion region for case E (10°, King) in Fig. [T8J We re- 
peat on the left the local polytropic index, while the right 
panel now gives the local vorticity value. Especially the lat- 
ter panel of Fig. [TB] hints at the presence of a sharp radial 
structure in the jet beam, itself surrounded by the instabil- 
ities of the mixing region. Fig. \W\ shows a radial cut of the 
jet at fixed axial distances of 2, 5 and 10 parsec for cases 
C and E (0°and 10°, King). We can see at the interface be- 
tween the jet and the shocked ISM that a layer forms with 
a width comparable with the beam width, with properties 
different from both media. For the cylindrical jet (opening 
angle of zero) this layer shows a lower density than the two 
surrounding regions with a clear local peak of the curl. For 
the 10 degrees opening angle we can see two such extrema 
in the local vorticity, present all along the jet from the root 
to the jet beam head, giving the beam a layered structure. 
These layered structures are sometimes invoked in studies 




Fig. 15. Details of two nodes along the axis. Internal struc- 
ture of the jet beam appears with a finite opening angle: a 
static shock, in the lab frame, forms along the path of the 
jet, which can be the seat of Fermi-type-II acceleration of 
particles. The size of the image is 7 x 8 pc with the two 
nodes at 3.5 and 7 pc along the symmetry axis. 
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Fig. 16. Comparison of the effective polytropic index for 
cases C (left) and E (right) at time 20 of the simulation. 
This figure reveals more clearly a series of secondary shocks 
that propagate in the shocked ISM region and also reheat 
it. 



of relativistic jets. lOstrowski fc Stawarzl (|200lD used this 
layered jet to explain increased emission at t he ex ternal 
border of the jet beam, whereas lO'Neill et~ al. (2005) used 
it as a way to keep the jet in equilibrium with its surround- 
ings. We recall that in our initial condition, the variation 
from the density of the jet to the density of the ISM is 
a step function. A sideways traveling rarefaction wave in- 
duces the shape we see eventually all along the entire jet 
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Fig. 17. Comparison of the effective polytropic index for 
cases D (left) and E (right) at time 20 of the simulation. 




Fig. 18. Presence of a dynamically formed layer surround- 
ing the beam of the jet. This is another zoom on the inlet 
(as in Fig.[T5]) for case E, now showing the (scaled) vorticity 
at the right. 



beam boundary at later times. Case E is more complicated 
because the beam is surrounded by the growing instabilities 
and because of re-collimation, the edge of the beam is not 
radially static in space anymore. 

In conclusion: a radially layered jet structure is present 
also without an opening angle for the injection of the jet, 
but the latter does amplify this structure in a finite rarefac- 
tion zone. 



5. Conclusion 

We quantified the fraction of energy present in the different 
jet-affected regions for realistic FR-II-type AGN jets. It was 
possible to show that the main mechanism of energy trans- 
fer from the jet to the ISM is through heating via shocks 
in the shocked ISM region, on top of mixing-related energy 
exchange. The energy transfer we measured for the SISM 
region dominates the same measurement for the mixing re- 
gion with a ratio of 3 to 1 on average. It is important to 
note that in reality, the mechanisms studied here coexist 
with magnetic energy transfers and radiation transfer. 

We discussed the role of finite opening angles on the 
dynamics of the jet, and in terms of efficiency of energy 
transfer. A wider opening angle increases the energy trans- 
fer but reduces the reach of the jet in its main direction 
of propagation and also reduces the portion of space that 
receives the energy. Our study followed relativistic Lorentz 
factor 10 jets for a lifetime of about 65 years and for a prop- 
agating distance up to 20 pc. We found that the amount of 
energy feedback from the jet to the ISM can vary from 10% 
to 70 % across all cases studied, where a total 10 46 erg.s" 1 
influx was assumed throughout. The jets were followed-up 
to asymptotic regimes where some of the energy transfer 
rates reached stationary values, which calls for larger scale, 
longer-time simulations to confirm these trends. This is be- 
cause we know that astrophysical AGN jets propagate on 
scales of kpc to Mpc for extended lifetimes. These longer 
time and larg e r scal e runs can be found in the works of 
iPerucho et al. (|201lh . They found that for early times, the 
energy transfer from the jet to its surroundings is domi- 
nated by heating by shocks before mixing processes become 
dominant for later times. This matches our results where 
for early times of the jet w e found a clear dom inance of 
the SISM heating by shocks lPerucho et all (|201lh adopted 
a resolution of 50 pc at 1 kpc from the source, assuming 
a cylindrical propagation of the jet making it complemen- 
tary to our own simulations. While their study reveals the 
long-term, long-distance trends with merely two gridpoints 
in the jet inlet beam, our study quantifies the effect of finite 
opening angles for early-time, beam-resolved cases. 

we also provided details on two dynamic phenomena 
due to finite opening angles. The first are the nodes cor- 
responding to recollimation shocks forming along the jet 
path. These nodes can be seats of significant acceleration 
of particles, which in turn sets a scenario to explain the X- 
ray nodes in observations of AGN jet s. This was originally 
studied in iKomissarov fe~F alle (1993), with the results of 
relativistic hydro simulations translated in synthetic im- 
ages. The opening angle of the jet plays an important role 
for the formation of these nodes, while the precise position 
and number of the nodes depends on the value of the open- 
ing angle. The second phenomenon we observed was the 
dynamically formed layered structure of the jet. Starting 
with a sharp profile for the density that separates the jet 
and its surrounding, the beam becomes surrounded by a 
clear density-structured layer of precisely the form needed 
to explain spine-sheath scenarios that are usually included 
for jet stability arguments or emission-related inferences. 
This feature does not need a finite opening angle to appear, 
but a higher opening angle does change its properties, such 
as a wider radial width and structured velocity profile. 

The next step of this work is to obtain synthetic X- 
ray maps of these simulations to compare them directly 
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Fig. 19. Radial cuts through the beam at distances z=2, 5 and 10 for an opening angle of 0° (left - case C) and 
10° (right - case E), both with a King atmosphere density profile. Top - density; center - curl of the velocities; bottom 
- axial velocity. The cuts reveal the position of the jet beam and of a layer between the jet beam and the shocked ISM. 
For the cylindrical case (left) the layer is near constant in the symmetry axis direction and has a width comparable to 
the radius of the beam. The opening angle case (right) shows a more complicated structure of this dynamically formed 
layer, with two peaks in the curl of the velocity. 



with observations. In this paper we intentionally focused on 
dynamical effects and energy transfer quantifications only. 
Moreover, we performed a parametric survey to identify 
clear trends due to opening angle variations and/or ISM 
density profiles. In future works, we intend to target specific 
well-diagnosed radio sources, and extend our study to the 
more computationally challenging 3D large-scale scenarios. 
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